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Abstract. We use Grassmann algebra to study the phase transition in the two- 
dimensional ferromagnetic Blume-Capel model from a fermionic point of view. This 
model presents a phase diagram with a second order critical line which becomes first 
order through a tricritical point, and was used to model the phase transition in specific 
magnetic materials and liquid mixtures of He 3 -He 4 . In particular, we are able to map 
the spin-1 system of the BC model onto an effective fermionic action from which we 
obtain the exact mass of the theory, the condition of vanishing mass defines the critical 
line. This effective action is actually an extension of the free fermion Ising action with 
an additional quartic interaction term. The effect of this term is merely to render the 
excitation spectrum of the fermions unstable at the tricritical point. The results are 
compared with recent numerical Monte-Carlo simulations. 
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1. Introduction 

The Blume-Capel (BC) model is a classical spin-1 model, originally introduced to 
study phase transitions in specific magnetic materials with a possible admixture of 
non- magnetic states [HE]. Its modification was also used to qualitatively explain the 
phase transition in a mixture of He 3 -He 4 adsorbed on a two-dimensional (2D) surface [3]. 
Below a concentration of 67% in He 3 , the mixture undergoes a so-called A transition: 
the two components separate through a first order phase transition and only He 4 is 
superfluid. On a 2D lattice representing an Helium film, He atoms are modelled by 
a spin-like variable, according to the following rule: an He 3 atom is associated to 
the value 0, whereas a He 4 is represented by a classical Ising spin taking the values 
±1. Within this framework, all the lattice sites are occupied either by an He 3 or 
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He 4 atom [3]. The 2D Blume-Capel model describes the behaviour of this ensemble 
of spins {S mn = 0, ±1}. In addition to the usual nearest-neighbour interaction, its 
energy includes the term A J2mn ^mw with S^n = 0, 1, to take into account a possible 
change in vacancies number. A can be thought as a chemical potential for vacancies, 
or as a parameter of crystal field in a magnetic interpretation. A simple analysis of 
the 2D BC Hamiltonian already shows that this model presents a rather complex phase 
diagram in the plane (T, A ), where T is the temperature in the canonical ensemble 
[I]. In the limit A — > — oo, the values S mn = are effectively excluded and the 
standard 2D Ising model is recovered, with its well-known second-order critical point 
at (T,A ) = (T c = 2/ln(l + V2) ~ 2. 269185, -oo), with the parameters taken in 
units of the Ising exchange energy J. At zero temperature, on the other hand, a 
simple energy argument shows that the ground state is the Ising like ordered state 
with \S mn \ = 1 if A < 2, and \S mn \ = else. There is therefore a first order phase 
transition at (T, A ) = (0, 2), suggesting a change in the order of the transition at some 
tricritical point at the critical line at finite temperatures. Mean field theory confirms this 
behaviour, and provides a second order transition line in the plane (T, A ) in the region 
extending from negative to moderated positive values of A [TJ El EJ IH E]. Beyond 
the tricritical point, as dilution increases, the transition becomes first order. Precise 
numerical simulations have been performed to study the phase diagram and to locate 
the tricritical point of the 2D Blume-Capel model [HI [TJ, El El ED]- From a theoretical 
background, several approximations have been used as well, such as mean field theory 
PQ El El HI E] , renormalization group analysis [HJ [12], and high temperature expansions 
[13J. Using correlation identities and Griffith's and Newman's inequalities, rigorous 
upper bounds for the critical temperature have been obtained by Braga et al [14J. It 
was also conjectured that exactly at the tricritical point the 2D BC model falls into 
the conformal field theory (CFT) scheme of classification of the critical theories in two 
dimensions [T5l [TBI ITT] . This is the case with m = 4 and c = 7/10, where c is the 
central charge [TBI [TTl [18] . The CFT analysis also implies a specific symmetry called 
supersymmetry in the 2D BC model at the tricritical point [TBI El EE] . 
The two-dimensional BC model is directly related as well to percolation theory [19] and 
dilute Potts model [20], where tricritical point properties are observed for percolating 
clusters of vacancies. We also mention quantitative results that match the universality 
class at the tricritical point of the BC model with the one of a 2D spin fluid model 
representing a magnetic gas-fluid coexistence transition [21], and similarities between 
BC phase diagram and Monte-Carlo results on the extended Hubbard model on a 
square lattice [22J. The advanced theoretical methods like bootstrap approach and 
perturbed conformal analysis, in combination with the integrable quantum field theory 
and numerical methods, have been applied to study the scaling region and the RG flows 
in the 2D BC universality class [281 124T 125] . 

The aim of this article is to present a different analytical method for the BC model 
in two dimensions with the use of the anticommuting Grassmann variables, originally 
proposed for the classical 2D Ising model in the case of free fermions [2~Bl |2~T] and 
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since then used to treat various problems around the 2D Ising model, such as finite 
size effects and boundary conditions [28], [29], quenched disorder [301 [31], an d boundary 
magnetic field [32j [33]. In contrast with the use of more traditional combinatorial and 
transfer-matrix considerations [3U [351 EH1 EH EH], this method is rather based on a 
direct introduction of Grassmann variables (fermionic fields) into the partition function 
Z in order to decouple the spin degrees of freedom in the local bond Boltzmann weights 
in Z. A purely fermionic integral for Z then follows by eliminating spin variables in 
the resulting mixed spin-fermion representation for Z. The method turns out to be 
particularly efficient to deal with models with nearest-neighbour interactions in the 
2D plane [261 EZ]- For the 2D Ising model, the fermionic integral for Z appears to 
be a Gaussian integral over Grassmann variables, with the quadratic fermionic form 
(typically called action) in the exponential (30J [38] . Respectively, the model is exactly 
solvable by Fourier transformation to the momentum space for Grassmann variables in 
the action. In physical language, this corresponds to the case of free fermions [381 [39] . 

As the additional crystal field term in the BC Hamiltonian is local, we hoped that 
the method will be applicable as well in this context. We will see in the following that 
though it is not possible to compute exactly the partition function and thermodynamics 
quantities of the BC model directly, since the resulting fermionic action for BC is 
non Gaussian, our approach allows to derive in a controlled way physical consequences 
from the underlying fermionic lattice field theory with interaction. In the continuum 
limit, a simplified effective quantum field theory can be constructed and analyzed in 
the low energy sector, leading to the exact equation for the critical line, that follows 
from the condition of vanishing mass, and to the effective interaction between fermions 
responsible for the existence of a tri critical point. The effects of interaction are assumed 
to be analyzed in the momentum-space representation. An approximate scheme such 
as Hartree-Fock-Bogoliubov (HFB) method can be used to locate the tri critical point. 
There are also some albeit formal analogies in this respect with the approaches typically 
used in the BCS theory of ordinary superconductivity. In general, it is interesting to note 
that in 2D a phase diagram of the BC model with first order transition and tricritical 
point can be described not only within a bosonic $ 6 Ginzburg-Landau theory [101 W\\ . 
where the order parameter is a simple scalar, but also with the use of fermionic variables. 

The article is organized as follows. After presenting the BC Hamiltonian and the 
related partition function in standard spin-1 interpretation, we apply the fermionization 
procedure leading to the exact fermionic action on the lattice. Then, from this result, 
we derive the effective action in the continuum limit and extract the exact mass. The 
condition of zero mass already gives the equation for the critical line in the (T, Ao) 
plane. The effective action also includes four-fermion interaction due to admixture of 
the S 2 = states (vacancies) in the system, with coupling constant g oc exp(— A), 
where —A = Aq/T, and Ao is the parameter of the crystal field in the Hamiltonian. 
We then give a physical interpretation for the existence of a tricritical point in the BC 
phase diagram by studying the fermionic stability of the BC spectrum at the critical 
line at order k 2 in momentum and compare our results with recent numerical Monte 
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Carlo simulations El U ED] • 

2. The 2D Blume-Capel model 

2.1. Hamiltonian and partition function 

The 2D BC model is defined, on a square lattice of linear size L, via the following 
Hamiltonian: 

L L L L 

H = — ^l'S'mn'S'm+ln + ^2>S'mn>Smn+l + A ^mn ■ (-Q 

m=l n=l m=l n=l 

In the above expression, S mn = 0, ±1 is the BC spin-1 variable associated with the mn 
lattice site, with m,n = 1,2,3, ...,L, where m,n are running in the horizontal and 
vertical directions, respectively. The total number of sites and spins on the lattice is L 2 , 
at final stages L 2 — > oo. The spins are interacting along the lattice bonds, Ji )2 > are 
the exchange energies. Notice that positive J\p > correspond to the ferromagnetic 
case. In addition to the Ising states with S mn = ±1, there are as well the non-magnetic 
atomic levels with S mn = 0, which we shall also refer to as vacancies. The crystal 
field parameter A plays the role of a chemical potential, being responsible for the level 
splitting between states S mn = and S mn = ±1. The Hamiltonian that appears in the 
Gibbs exponential may be written in the form: 



L L L L 

— (3H = ^iS' mn S' m+ i n + K2S mn S mn+ i + A , (2) 

m=l n=l m=l n=l 

where = (3Ji,2 are now the temperature dependent coupling parameters, (3 — 1/T 
is inverse temperature in the energy units, and A = — (3A . In what follows we will 
assume, in general, the ferromagnetic case, with positive Ji_2 > and > 0, though 
the fermionization procedure by itself is valid irrespective of the signs of interactions. |ij 
The positive A (negative Ao) is favourable for the appearance of the Ising states with 
S^ n = 1 in the system, with the ordered phase below the critical line in the (T, A ) 
plane, at low temperatures, while negative A (positive A ) will suppress Ising states, 
being favourable for vacancies. In the limit A — > oo, or Ao — > — oo, the states with 
^mri = ® are effectively suppressed and the model reduces to the 2D Ising model, with 
the critical temperature being defined by the condition sinh 2K\ sinh 2K 2 = 1. As Ao 
increases to finite values, there will be a line of phase transitions in the (T, A ) plane. 
The increasing A admits the appearance of the vacancy S^ n = states. Respectively, 
the critical line goes lower as Ao increases from negative to positive values and terminates 
at A = J\ + J 2 at zero temperature, so that all sites are empty at larger positive values 

\ In what follows, by presenting the numerical results, we shall typically assume isotropic case for 
interactions in the above Hamiltonians, with J\ = J2 = J, K\ = K2 = K, and K — (3 J. We will also 
use, in some cases, the dimensionless parameters normalized by the exchange energy J for temperature 
T and chemical potential Aq. 
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of A at T = 0. A remarkable feature of the BC model is that there is also a tricritical 
point on the critical line at finite temperatures somewhere slightly to the left from 
Ao = J i + J2, where the transition changes from second to first order. 
The partition function Z of the BC model is obtained by summing over all possible spin 
configurations provided by {S mn = 0, ±1} at each site, Z = J2s=o ±1 e_/3H — ^ r {s} er^ H ■ 
Using the property {S mn = 0, ±1}, it is easy to develop each Boltzmann factor appearing 
in the above trace formula in a polynomial form: 

exp(KiSS') = 1 + XiSS' + X'^S 12 , i = 1,2, (3) 

with 

Xi — smhKi, X' i = coshK i — l. i = l,2. (4) 

The partition function is then given by the product of the above spin-polynomial 
Boltzmann weights under the averaging: 



L L 



Z= Tr 

{5 m „=0,±l} 



ACJ2 



{nn 

m=l n=l 

X (1 + A 2 S mn S mn+ i + X' 2 S'mn^mn+l) j 



v,l "I - -^i S mn S m +i n + XiS mn S m ,i 



(5) 



This expression will be the starting point of the fermionization procedure for Z using 
Grassmann variables we develop in Section 3. At first stage, we introduce new 
Grassmann variables to decouple the spins in the local polynomial factors of expression 
(J5j). At next stage, we sum over spin states in the resulting mixed spin-fermion 
representation for Z to obtain a purely fermionic theory for Z. 



2.2. Local spin decomposition 

In what follows, we shall need to average partially fermionized Z over the spin states at 
each site. This averaging will be performed in two steps, first we keep in mind to average 
over the Ising degrees of freedom, S mn = ±1, then adding the contribution of vacancies, 
Smn — 0- The two cases may be also distinguished in terms of variable 5*^ n = 0, 1. In 
this subsection, we shortly comment on the formalization of this two-step averaging. 
Provided we have any function of the BC spin-1 variable f(S mn ), with S mn = 0, ±1, the 
averaging rule is simple: 

f(Smn) = f(0) + f( + l) + f(-l). (6) 

Smn — 0, it 1 

In forthcoming procedures, we ought to average first over the states S mn = ±1 at each 
site, provided S^ n = 1, while making the sum over choices S^ n = 0, 1 at next stage. In 
principle, since S mn = sign{5' TOn } \S mn \, with sign{S mn } = ±1 and \S mn \ = S 2 m = 0, 1, 
we can try simply to write S mn = y mn P mn , where y mn = 0, 1, and a mn = ±1, and to 
average over the component states y mn = 0, 1 and a mn = ±1 as independent variables. 
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This gives: 

fiVmn^mn) = /(+0) + /(-O) + /(+1) + /(-l) . (7) 

ymn —0,1; CTm.n — dz 1 

We see that the zero state is counted twice, in contradiction to (jSJ). This may be 
corrected by introducing in the definition of the averaging the weight factor | at y m „ = 0. 
Equivalently, this may be formalized by adding 2~ 1+Vmn under the sum. This results 
the sum of three terms in agreement with (jHJ): 

J- 2- 1+ y- f(a mn y mn ) = f(0) + /(+1) + /(-l) . (8) 

In fact, this decomposition scheme with S mn = a mn y mn and independently varying 
a mn = ±1 and y mn = 0, 1 is somewhat more close to the situation for the two- 
dimensional Ising model with quenched site dilution [30], [31]. In that case a mn = ±1 is 
simply the Ising spin, while the variable y mn = 0, 1 is the quenched dilution parameter, 
counting whether the given site is occupied or dilute, and both averaging rules (JTj) and (JH]) 
can be interpreted physically. The case (JJ|) means in fact that there is a spin <j mn = ±1 
also at site y mn = 0, which is not interacting with its nearest-neighbors. This empty, or 
rather disconnected site, by flipping over two states ±1 under temperature fluctuations, 
will give however a contribution to the entropy, In 2 by empty site. The case (jHJ) means 
that the site y mn = is really dilute, or empty, with no spin degree of freedom at it, 
even disconnected. For the quenched dilute 2D Ising model, the quenched averaging 
over some fixed temperature-independent distribution y mn = 0, 1 is physically distinct 
from the a mn = ±1 averaging, and is assumed to be performed rather on — f3F = In Z, 
but not on Z itself. The situation is different for the BC model, which is in essence the 
annealed case of the site dilute Ising model, with the averaging simultaneously over all 
states S mn = 0, ±1 at each site for Z itself. In this case the averaging is to be performed 
strictly according to the rules like ([6]) and (jHJ), but not ([7j). 

There is still another way to formalize the averaging over the possibilities of 
S mn = ±1 before we actually fix S^ n = 0,1. It is based on the observation that 
the result of the averaging (jSJ) will not be changed if we replace S mn — > a mn S mn , with 
a mn = ±1, since the sum includes S mn = ±1 anyhow: 

f(S mn )= Yl f(°mnS m n) = f(0) + f( + l) + f(-l), O mn = ±1 . (9) 

Sjnn =0,zbl Srnn = 0,il 

Thought the above equation holds already for any fixed value of a mn = ±1, we can as 
well average it over the states a mn = ±1, introducing factor | for normalization. The 
averaging of f{a mn S mn ) itself gives: 

2 f{ a rnnS mn ) = -\f{S m n) + f{-Srn n )}= g(S 2 mn ) , ^ = 0,1. (10) 

&mn — il 

The result of the averaging will be a function g which only depends on IS^I = 0, 1, alias 
S^ n = 0, 1, but not on the sign{S' mn } of S mn = sign{5' mn }|5' mn |. In terms of g(S^ n ), 
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the equation results: 

E f( S mn)= {\ E /M„)}=5(0) + 2 9(1). (11) 

Smn — 0, i 1 Srnn — 0, i 1 (Tmn— il 

In this form the two-step averaging will be realized in the procedure of elimination of 
spin variables by constructing the fermionic integral for Z in the forthcoming discussion. 

3. Fermionization and lattice fermionic field theory 

The expression of the BC partition function Z as a product of spin polynomials under 
the averaging as given in ([5]) will be the starting point of the fermionization procedure 
for Z. This procedure has first been introduced in the context of the 2D pure Ising 
model [261 EZj- It relies on interpreting each spin polynomial Boltzmann weight in ([5]) 
as the result of integration over a set of two Grassmann variables, which decouples the 
spins under the integral. Before going into details, we remind in the following subsection 
few essential features about Grassmann variables and the rules of integration. 

3.1. Grassmann variables 

Mathematically, Grassmann variables may be viewed as formal purely anticommuting 
fermionic numbers [12]. In physical aspect, they are images of quantum fermions in 
path integral [12]. We remind here few basic features about Grassmann variables that 
are needed in the rest of the paper. More details can be found in [131 44j . A Grassmann 
algebra A of size N is generated by a set of N anti-commuting objects {a\, a 2 , . . . , a^} 
satisfying: 

OjOj + OjOj = , a\ = , i, j = 1,2, . . . , N . (12) 

This as well implies including the case i = j. Unlike quantum fermions, 

Grassmann variables are totally anticommuting. Note that any linear superpositions of 
the original variables (|T2|) are again purely anti-commuting with each other and with the 
original variables, and their squares are zeroes. Functions defined on such an algebra are 
particularly simple, they are always polynomials with a finite degree (since of = 0). It is 
possible to define the notion of integration [1TI W2\ H3l HI] in algebra of such polynomials 
with the following rules. For one variable, the rules are: 

J daj • aj = 1 , J daj -1 = 0. (13) 

The integral with many variables is considered as a multilinear functional with respect 
to each of the variables involved into integration (integral of a sum is the sum of the 
integrals etc). In multiple integrals, the fermionic differentials are assumed again anti- 
commuting with each other and with the variables themselves. The integration of any 
polynomial function of Grassmann variables like f(a) — f(ai, a,2, ■ ■ ■ , oat) then reduces, 
in principle, to a repeating use of the above rules. The rules of change of variables 
in Grassmann variable (fermionic) integrals under a linear substitution are similar to 
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the analogous rules of common (commuting) analysis. The only difference is that the 
Jacobian of the transformation will enter now in the inverse power, as compared to the 
commuting (bosonic) case [HJ 1121 B21 S3] • 

With the above definitions, the Gaussian integrals over Grassmann variables are all 
expressed by the equations relating them to determinants and Pfaffians. The basic 
equation for the determinantal integral of first kind reads: 



where the integration is over the doubled set of totally anti-commuting variables {a, a*}. 
The (square) matrix A in the exponential is arbitrary. In applications, the quadratic 
fermionic form in the exponential like in (fl4"l) is typically called action. Since the 
action is quadratic, the integral is Gaussian. The exponential in (fT4l) is assumed in the 
sense of its series expansion. Due to nilpotent properties of fermions, the exponential 
series definitely terminates at some stage, thus resulting a finite polynomial in variables 
involved under the integral. [With respect to the action S = aAa* taken as a whole, 
the last nonzero term will be with S N ^ 0, while S N+1 = 0. Alternatively, the same 
polynomial for the exponential from ( f!4l will follow by multiplying elementary factors 
like exp(aiAija*) = 1 + a,iAija*}. In physical interpretation, the integral of the first kind 
( TT4|) with complex-conjugate fields rather corresponds to Dirac theories. 
The Majorana theories with real fermionic fields are presented by the Gaussian integrals 
of the second kind related to the Pfaffian. The basic identity for the fermionic integral 
of the second kind reads: 



The integration is over the set of even number N of Grassmann variables, the arrow 
in the measure indicates the direction of ordering of anti-commuting differentials. The 
matrix in the exponential is now assumed skew-symmetric, Aij+Aji = 0, An = 0, which 
property is complimentary to fermionic anticommutativity. The result of the integration 
is the Pfaffian associated with the skew-symmetric matrix A from the exponential, 
otherwise, one can associate the Pfaffian on the r. h. side of (TT5"j) with the above-diagonal 
triangular array of elements of that matrix, {Aij | 1 < i < j < N}. In mathematics, the 
Pfaffian is known as a certain skew-symmetric polynomial in elements of a triangular 
array of the above kind. In physics, the combinatorics of the Pfaffian also is known 
under the name of the (fermionic) Wick's theorem. Note that the identity (fT5|) can be 
assumed by itself for the definition of the Pfaffian. 

In a combinatorial sense, the determinant is rather a particular case of the Pfaffian. 
Respectively, the integral ( fT4l) is a subcase of the integral ( fTBT) . It can be shown, on the 
other hand, that (PfA) 2 = det A for any skew-symmetric matrix A. This implies that, 
in principle, an integral of the second kind (fl5l) can always be reduced to an integral 
of first kind (fl4l) by doubling the number of fermions in ( 1151) . In applications like in 




(14) 




(15) 
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the Ising and BC models, where the original integrals in the real lattice space rather 
appear in the Pfaffian like form of (j!5p . this reduction to the determinantal case occurs 
automatically after the transformation to the momentum space, where the fermionic 
variables are typically combined into groups of variables with opposite momenta (k, —k), 
which play the role of the conjugated variables like in (114p . In practice, for low- 
dimensional integrals, most of calculations can be performed simply from the definition 
of the integral, by expanding the integrand functions into polynomials. 

3.2. Fermionization procedures 

In the same spirit as for the 2D Ising model [27], we introduce two pairs of Grassmann 
variables (a mn ,a mn ) and (b mn ,b mn ) to factorize the polynomials appearing in (JSJ). 
Namely we use the relations: 

i + VQ S , i 4-A'^ 2 S 2 — do dn p^+K^nS^J a mn a mn 
1 f ^■l>Jmn' J m+ln > A l> J mn lJ m+ln — / uu mn uu mn e 

X (1 + Ct mn S mn ) (1 + Ai a mn >S'm+ln)) 

1 4- \oS S , 1 4- A\S" 2 .S" 2 — db dh pC+^^^L-fi)^™^" 

x(l + 6 mil Smn) (1 + ^2 &mnSmn+l)- (16) 

For the sake of simplicity in notation we introduce the following link factors: 

We also define the Grassmann local trace operators which associate to any function 
/(...) on the Grassmann algebra as follows: 

n o(WAX+l«) a ranVTl f/n 

mn; u 'mnjj 



2<t [/("rani ^mn)] / d(Z mr jdG m?1 1 mn m+1 "^ mn mn /(fl mn ,(3 

It [ffft ft )] = dh db p( 1 + A 2 s mnS'^, l+1 )f) m „S mn --/i t \ /-,^ 

\ u mni u mn)\ — / { - xu rnn'~ llJ inn c J \ u mni u mn) ■ \ LO ) 

(bmn) J 

The factorized Boltzmann weights from (fT6l) now read: 



1 + A 2 5'^ ri S'^ n+1 + A 2 5'^ n 5'^ ri+1 = Tr [S mn 5 mn+1 ]. (19) 

(^777,77, j 

Introducing the above Boltzmann weights into the original expression (jSJ) for Z, we 
obtain a mixed representation containing both spins and Grassmann variables for Z. 
Notice that as the separable link factors like A mn , A mn , B mn , B mn are neither commuting 
nor anti-commuting with each other, the order in which they appear in the product may 
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be important. The factorized bond weights, however, presented in (1191) by doubled link 
factors under the trace operators, are totally commuting, if taken as a whole, with any 
element of the algebra under the averaging. For the whole lattice, following the rules 
([18]) . we define the global trace operator as follows: 

» L L 

5t- \f] / I I I I d < ^mndcimnd& ?T md& mrl C m ™ j ' iflmni O'mni^mni^mn) (^0) 
(a,6) J 



in 1 n=l 

L L 



X 



ex P I Z-^ ^l^mn^m+ln) a mii a rrm + (1 + ^2^71^)7111+1)^"^] f ■ 

l.m=l ra=l J 

The all even-power terms in spin variables are now incorporated into the generalized 
Gaussian averaging measure of (1201) . including the term with chemical potential. The 
partition function is then given by 



Z = Tr It 

{S} (a,b) 



n n 

n=l \ m=l 



i an-^-m+ln) (-Smn-^mn+l)) 



(21) 



At this stage the factorized partition function appears as a double trace, over the spin 
degrees of freedom, with Tr{s}, and over the Grassmann variables, with Xt( aj /j) . The 
idea of the next step is to make spin summation in ()2ip to obtain a purely fermionic 
integral for Z. At first stage, we keep in mind to eliminate rather the Ising degrees ±1 
of spin variables in (l2~Tj) . The averaging over S^n = 0, 1 will be performed at next stage. 



3.3. The ordering of factors 

Up to now we only add extra fermionic (Grassmann) variables to obtain the mixed 
expression (j2ip . where the spin variables are actually decoupled into separable link 
factors like ( jTTl) . Further algebraic manipulations are necessary to simplify this 
expression in order the spin averaging be possible in each group of factors with the 
same spin. For any given mn, there are four such factors, A mn , B mn , A mn , B mn , which 
all include the same BC spin S mn = 0, ±1. What we need is to be able to keep nearby 
the above four factors with the same spin at least at the moment of the spin averaging. 
We apply the mirror-ordering procedure, introduced originally for the two-dimensional 
Ising model, to move together, whenever possible, the different link factors containing 
the same spin. Despite of that the separable link factors like (1171) are in general neither 
commuting nor anticommuting, it is still possible to make use of the property that the 
doubled combinations like A mn A m+ i n and B mn B mn+ i are effectively commuting, if taken 
as a whole, with any element of the algebra under the sign of the Gaussian fermionic 
averaging in (I2TI) . Using the notation for the ordered products similar to that of [26| [2T] . 
this leads to: 

L L 

Z = Tr Tt | J^J J^J ^(A mn A m+ i n )(B mn B mn+ i) \ 

' m=l n=l 
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L L l, 

Tr Tr •{II II B mn A mn A m+ i n ■ I I B r 
{S} (a,b) I A A L A A A A 



n=l m=l m=l 



Tr £r 

{S} ( a ,6) 




i-^mn^mn I ' I J^J B, n 



ran 1 

n=l \ \ m=l / \ m=l 



(22) 



In the above transformations, we use mirror-ordering decoupling for factors in vertical 
direction, B mn B mn+ i, with respect to n, then insert the commuting factorized 
horizontal weights, A mn A m+ i n , and reread the resulting products in few subsequent 
transformations (cf. [26j [2?]). The boundary terms are also to be treated properly as 
we pass from (0) to ( 1221) . The simplest case is provided by the free boundary conditions. 
The free boundary conditions for spin variables, Sl+iu = S m L+i = 0, in (jSJ) correspond 
to the free boundary conditions for fermions, a 0n = b m0 = 0, in ( |22l) . For free boundary 
conditions, the transformation from (JHJ) to ( 1221) is exact. In what follows, however, we 
will typically assume the periodic boundary conditions for fermions in representations 
like (1221) . These are most suitable closing conditions when passing to the Fourier space 
for anticommuting (Grassmann) fields. The change of the boundary conditions of this 
kind is inessential in the limit of infinite lattice as L 2 — > oo. In principle, one can pay 
more attention to the effects of the boundary terms in the periodic case, which can 
actually be treated rigorously also for finite lattices [271 1211 [291 [32j [33] . 
In the case of the 2D Ising model, with S^ n = 1, we can explicitly perform the trace over 
the Ising spin degrees of freedom S mn = ±1 recursively at the junction of two m-ordered 
products in the final line of ( 1221) . The situation is slightly different in the BC case, since 
S m n = 0) I; instead. Also, the trace operator (|20|) contains terms with S^n = 0, 1 which 
are coupled at neighboring sites. Therefore it is not possible to trace over the whole set 
of states S mn = 0, ±1 in the BC case directly in ( 1221) . but only we can eliminate first 
the Ising degrees sign{S mn } = ±1. The BC variables S 2 mn = 0, 1 will still remain as 
parameters and will be eliminated at next stages. The elimination of the Ising degrees 
we will realize by the symmetrization transformation, S mn — > <J mn S mn , with averaging 
over a mn = ±1, following the procedures explained in (l9l)- (|TTl) above. The details of the 
a mn = ±1 averaging are discussed in the next subsection. Thus, the ordering procedure 
on the link variables (TlTl) allows us to eliminate at least a one part of the spin degrees 
of freedom in the factorized expression for Z resulting in the final line of ( |22l) . 



3.4- Spin summation 

At the junction of the two ordered products in ( |22l) . with S mn — > cr mn S mn , we perform 
the trace a mn = ±1 recursively, for m = L, L—l, . . . , 2, 1, for given fixed n, starting with 
m = L. The procedure will be then repeated for other values of n — 1, 2, 3, . . . , L. The 
four relevant factors A mn , B mn , A mn , B mn with the same spin that met at the junction of 
the two m-product in (1221) . for given n, are to be specified from (jXTJ) . There we assume 
S m n o~ mn S mn . Then we multiply the above four factors, taking into account that 
a mn = 1, so that S^n ~^ a mnS m n ~^ ^mn' an d sum over the states a mn = ±1. This will 
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eliminate all odd terms in the polynomial so obtained. The averaging thus results: 
- \ A R 4 R 

2 / r ^mn- LJ mn- n -mn- LJ mn 

1 "I - S mn CL mn b mn -\- S mn (\i CL m —ln A2 bmn—l)(flmn ^mn) "I - *S'rrm / ^l / ^2 
^rnn^^^2&mnbmn&m—lrfirnn—l 
= exp S^j„ (a mn b mn + (Ai(3 m _ ln + A2& m n-l)(°W + &mn) + Ai^Om-ln^rTm-l^ • (23) 

The even fermionic polynomial resulting under the averaging can be written as a 
Gaussian exponential, as is shown in the final line. This term is totally commuting 
with all other elements of the algebra and can be removed outside from the junction. 
The BC spins still remain in the form of S^ n = 0, 1 in (I2"3"j) . but the Ising degrees, 
sign{S mn } = ±1, are already effectively eliminated. After completing the above 
averaging procedure at the junction at m = L, for given n, we repeat the calculation for 
m = L — 1, . . . , 2, 1, and then for other values of n = 1, 2, . . . , L. Adding the diagonal 
terms from the definition of the fermionic averaging (J20l) . the partially traced partition 
function finally reads: 

f L 

Z = 2 TV / da mn da mn <lb mn db 

mn exp [AS'j^yj + (1 + A 1 S mn S m+ i n )a mn a mn 

{5 2 =0,1} J " 



m,n=l 

->2 



+ (1 + A 2 S mn S mn+l )b mn b mn + S mn (Aia m _i„ 

mn—1 1 \(J"mn i "mn ) 

^mn dmnbmn S mn \\ \2Clm— Inbmn— l] • (24) 

The resulting integral for Z in ( 1241) is the Gaussian integral, which includes yet the 
variables S^ n = 0, 1 as parameters. At this stage, it easy to recognize that the 2D Ising 
model is solvable, since in this case S^n = 1 at all sites. The partition function Z is 
then given by a Gaussian fermionic integral, which can be readily evaluated by passing 
to the momentum space for fermions [26j [27j . This results the Onsager's expressions 
for Z and —f3F = In Z. In the BC model case, it remains yet to eliminate S"^ n = 0,1 
degrees of freedom in the above expression ([21]) for Z. 

The trace over S^ in = 0, 1 can be performed in (12^1) after we manage to decouple the 
variables in terms including S^ nn S^ n+ln and S'^ n S'^ n+1 . Several methods are possible. 
A one way is to introduce another auxiliary set of Grassmann link variables, similarly 
to what we previously did to decouple the factors S mn S m+ i n and S mn S mn+ i in QT5|) . 
It is possible however to avoid the introduction of the new fields by using instead the 
following rescaling of the fermionic variables under the integral: a mn — > a mn / S^„, b mn — > 
bmn/Smn- Respectively, to preserve the integral invariant, one has to rescale the 
differentials: da mn — > S^da^, db mn — > S^ n (i6 mn . This may be viewed, in principle, 
as a kind of a change of variables in a fermionic integral, which leaves the integral 
invariant, as it follows from the basic rules of integration. The variable S^ nn then 
disappears in some places inside the exponential and appears in the others, the terms 
with S^ n S^ +ln and 5'^ n S'^ n+1 being decoupled. Also, the resulting seemingly singular 
expressions like S^ n exp(a mn a mn / S^) are to be understood as exp(a mn a mn / S^ n ) — 
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^(i + V^mn/^J = S^ n + a mn a mn . Finally, after shifting some indices in the sums, 
we obtain: 

f L 

Z 2 Tr / II do mn dci mn d& ?Tm d6 mn (>S , jrm -\- QjmnQ>mn)(S mn ~\~ b mn b mn ^ 

{S 2 =0,1} J 

m,n=l 

X exp [AS^ + S'^ nn [\[ a m -\n ^m-ln + A 2 b mn -lbmn~l + A1A2 flm-ln^mn-l)] 

X exp [a mn b mn + (Aia m _i„ + \2bmn-l) {&mn + &mn)] • (25) 

In this expression, we can already locally perform the sum over S^ n = 0, 1 at each 
site. The rules like ()8l)- (jTT!) are to be taken into account in order not to count twice 
the contribution of S^ nn = states. By averaging the part of the product explicitly 
depending on S^n = 0, 1, we obtain: 

{5^=0,1} 

^m—ln^m—ln 

+ Kb 

mn—lbmn—l Q"m—lnbmn—l ) f 

Q , mnQ"mnbmnbmn 26 6 , (26) 

where G mn in the exponential in final line stands for the local part of the action resulting 
at the Ising site with S^ nn = 1. The first term in final line is the one produced at dilute 
site with S^ n = 0. The explicit expression for G mn reads: 

G m n Q"mnQ"mn b mn b mn -\- A1A2 Q"m—lnbmn—l A^ d m —.\ n d m —\ n ~\~ A2 b mn — \b mn — \. (27) 

The result of the averaging from (1261) can as well be written as a unique exponential 
taking into account the nilpotent property of fermions: 

a a k k -4-9 p A p Gmn — 9p^p Gmn ( 1 -I- — a a h b p — A-G mn 

Uj mn Uj mn u mn u mn ~ * io \ 2 n n 

26 exp (^G mn -\- ~d mn tt mn b mn b mn 6 j 

2c exp ^G rnn -\- —& &mn&mnbmnbrnn& mn ^ . (28) 

In the final expression, we assume the local action G mn to be replaced by its reduced 
version G' mn , since the prefactor a mn a mn b mn b mn annihilates the two first terms of G mn 
in the exponential. The reduced action reads: 

+ A1A2 (29) 

Substituting this result into (125]) and shifting the mn index in some of the diagonal 
terms of the resulting combined action, we obtain: 

/L L L L 

\\ \\ da mn da mn db mn db mn exp j ^ ^ [ (1 + A'J a mn a mn 
m=l n=l m=l n=l 

^2) bmnbmn ^mnbmn (AlO-m— In ^2^mTi-l)(^mn b mn ) 

-I-A1A2 90 ^"mn^mnbmnbmn exp ( A^ Qj m —\ n & m —\ n A2 b rnn —\b mn —\ 



— A1A2 a 



m—lnUmn—1 1 



9o = e- A /2, (30) 
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with —A = +/3A . This is already a purely fermionic integral for Z, the spin degrees of 
freedom being completely eliminated. The interaction part of the action is introduced 
with coupling constant go oc e /3A °, which depends on Ao- To simplify the comparison 
with the 2D Ising model, and for other needs, we now rescale some of the Grassmann 
variables under the integral using the following transformation: 

(1 + Aj) tt mn > d-mri} (1 4" A 2 ) b mn > b mn . (31) 

The corresponding differentials are to be rescaled with inverse factors. In this way, we 
obtain the final result of this subsection: 

/L L L L 

II II da mn da mn db mn db mn exp j ^ ^ 
m=l n=l m=l n=l 

QimnQ'mn 4~ b Tnn b mn 4" CL mn b mn 4~ (t\d m —\ n 4" t'2^mn— 1 ) ("run "I - b mn j 4" tl^2 A'm-ln^mn-1 
"I - fl'O ^mn^mnbmnbmn GXp ( 'y^Cl m _ \ n d m — \ n 72^mn-l^mn-l ^1^2 ^m-ln^mn-l) } > (32) 

where £j = tanh ifj and we have introduced the following constants: 



-A 



1 



^-2coshir lC oshir 2 ' 7i - 1 "^^~ 1 ~V 1 -^ (33) 
In a compact form, the equation (|32j) reads: 

Z = (2e A cosh if i cosh iv~ 2 ) L2 ^ DaDaDbDb exp(S Ising + <S int ) . (34) 

Note that the fermionic integrals (130]) and (|32|) for the Blume-Capel partition function 
Z are still the exact expressions. They both are equivalent to each other and to (jSJ). The 
above correspondence is exact even for finite lattices, provided we assume free boundary 
conditions both for spins and fermions. [jj We can recognize in (132]) and (1341) the Ising 
action, which is simply the Gaussian part of the total action (cf. [261 I2T]): 

L 

Rising ^ ^ Q"mnQ"mn 4~ bmnbmn 4~ ®"mnbmn 



m,n=l 



+ (tl )(a 

inn 4- b mn ) + tit 2 a m 

— lnbmn—li 

(35) 

and the non-Gaussian interaction part of the total action, which is a polynomial of 
degree 8 in Grassmann variables (after expanding the exponential): 



C _ n \ n n h h p —H a m-\na. m -ln—12ilmn-lbrnn-l—t\t2a m -lnb mrl -l ( "if{\ 

t -'int !J0 / Uj mn Uj mn u mn u mn' : -' • 



L 

E 

m,n=l 

The BC model differs from the Ising model by the interaction term ( 136]) in the total 
action, which is not quadratic. Therefore the BC model is not solvable in the sense of 
free fermions, as distinct from the pure 2D Ising model. 

§ The exact fermionic representation for Z for a finite lattice with periodic boundary conditions for 
spin variables in both directions also can be derived. The result will be the sum of four fermionic 
integrals like (|30|) and ()32|) . with periodic-aperiodic boundary conditions for fermions, in analogy to 
the case of the 2D Ising model on a torus (cf. [27l [28j [29] ) . 
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It may be still of interest to try to recognize the structure of the phase diagram of the 
BC model directly from the fermionic integrals fl30|) - fl36|) before actual calculation. The 
interaction is introduced in the above BC integral ( 132]) for Z with the coupling constant 
go oc exp[/3(A — J\ — J2)]. In the limit A — > 00 (or A — > —00), which corresponds to 
go = 0, the gap between the two degenerate states S — ±1 and the singlet state S = 
becomes infinitely large and the model reduces effectively to the 2D Ising model. For A 
finite, the coupling constant go is finite and the presence of the vacancy states becomes 
possible. The coupling constant go increases as the number of the vacancies in a typical 
configuration of a system increases, with increasing A . At zero temperature, on the 
other hand, as (3 — > +00, we find g = for A < J\ + J2, which corresponds, again, 
to Ising ground state, while for A > J\ + J2 we have go — > +00, which will mean that 
the ground sate is empty (all sites are occupied by vacancies). These features at T = 
can be readily guessed already from the form of the original Hamiltonian (DO). A more 
sophisticated analysis of the integral (1321 will be needed to define the precise form of 
the critical line and to locate the tricritical point at that line with increasing dilution. 
In the following, for simplification, we will only consider the isotropic coupling case, 
with K\ = K 2 = K, ti = £2 = t and 71 = 72 = 7. 

3.5. Partial bosonization 

The previous action contains two pairs of Grassmann variables per site. This can not 
be reduced to a one pair (minimal action) unlike the Ising model, where half of the 
variables are irrelevant in the sense they do not contribute to the critical behaviour and 
can be integrated out already at lattice level. The point is that the reduced Ising action 
with two variables per site readily admits QFT interpretation and simplifies the analysis 
in the momentum space [301 SB ES [39]. In the BC case, the two pairs of fermions are 
coupled together by Eq. (13"6i) . preventing a direct integration over extra variables like 
a mn , b mn . However, as we will see in the following, it is still possible to recover the 
minimal Ising like action with a one pair of fermions per site using auxiliary bosonic 
variables. In the interaction part of the action (I3"6j) . it is indeed tempting to replace the 
products a mn a mn and b mn b mn , which are formally looking similar to occupation number 
operators, or local densities, by the new commuting variables as follows: 



Vmn Q'mnQ'mn > 7~mn U mn U mn , Vmn ^~mn ^ • \"'J 

These new variables rj mn , r mn are nilpotent (as Grassmann variables) but purely 
commuting: that is why we will abusively call them (hard core) "bosons". In the 
following, we will add also one more pair of commuting nilpotent fields fj mn ,f mn , to 
put integrals into a more symmetric form. The identities like ( 1371) are rather to be 
understood in the sense of correspondence, to be realized by a properly introduced 
delta-functions (Dirac distributions). This eventually allows us to reduce the degree of 
polynomials in Grassmann variables by a factor 2 each time the replacement like (I37|) is 
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performed, even if terms like a m _i n 6 mn _i in (}3"oT) can not be replaced. We will see below 
that, in principle, we can write down an action containing a one pair of Grassmann 
variables and a one pair of bosonic ones per site. To do so, we introduce the following 
Dirac distribution for any polynomial function / of nilpotent variables like a mn a mn or 



f(a mn a mn ) - J d Vmn dfj mn f( Vmn ) exp [fj mn (Vmn + a mn a mn )} , 

We assume a natural definition of the integral for commuting nilpotent variables with 
the following rules of integration (similar rules are assumed for fj mn , f mn ): 

drjmn (1, Vmn) = (0, 1) , / dr mn (1, T mn ) = (0, 1) . (39) 



For application of the rules like ( 1391) in the QFT context also see [35]. Applying ( |38l ) 
directly in (I32I) . we obtain the integral with the following expression for the action: 

m,n 

~^~T)mn\Vmn ^mn^mn) 7~ mn {T mn ~\~ b mn b mn ) . (40) 

We can now integrate over the a mn 's and 6 mn 's, and replace formally, for convenience, 
the variables a mn by c mn and b mn by — c mn in the remaining integral. We obtain: 

L 

<S ^ ^ ^C m?1 C mn (l -\- T mrl )(l -|- Tj mn ) -\- Tj mn Tj mn -\- T mn T mn 



mn=l 



_ l _ [Cmn(l Vmn) C mn (l -|- T mn )]t(c m _ \ n -\- C mn — i) t C m —\ n C mn —\ 



mn 



Ifijlm—ln ^mn-l) "I - T Vm—ln^mn—l t ^m—ln^mn—1 f • (^1) 



The advantage is that now there are only two fermionic variables per site, which is 
suitable for the QFT interpretation [301 EH]- Note that the integral associated with 
the action ( HT1) will still be the exact expression for Z. The number of the fermionic 
variables being reduced, the next operation is to try to integrate out, whenever possible, 
the auxiliary bosonic fields from action (141|) . In fact, we can further integrate over one 
pair of bosonic variables, for example r mn , f mn , using the integration rules like fl38|) . 
since 



dT rnn d,T rnn j '{Tmn) 6Xp[T mn (T Jnn t(c m —\ n C mn —\}C mn ~\~ C mn C mn (l -|- ^mn))] 

f [ ^(Cm— In C mn „x)c mn -|- C mn C rnn {\ -\- T] mn )~\ ■ (^^) 

There f{T mn ) may be any function of nilpotent variable r mn . We could also have chosen 
to integrate over the r) mn ,fj mn instead. Integrating over T mn ,f mn according to ( 1421) . we 
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finally obtain the reduced integral with the local action: 

& CmnCmn ^(^'rmi ^mn){pm—ln Cmn—l) t ^m—ln^mn—1 
~^~T}mnT1mn T}mn Cmn ^(c»n- In ^mn-1, 

m—ln + Q mn—1 ) 

~^~ r ) Vm—lnQmn—1 t C m —\ n C mn —\ , (43) 

with 

Qmn [Cmn(l "I - T)mn) ^(Cm- In Cmn— l)]c»nn . (44) 

It is easy to recognize in the first line of ( 1431) the minimal local action for the pure Ising 
model (30J, EH] with one pair of Grassmann variables per site: 

Rising C mn C mn -\- t(c mn -|- C mrl )(c m _i n C mn — i) t C m _i n C mrJ _ l- (45) 



This is the same action that follows by integrating a mn , b mn from ( 1351) . The rest of the 
action describes the interaction between fermions and bosons: 



^int TjmnQmn f]mn^n 



Cmn t{fim—ln Cmn— 1, 



+ go VmnQmn 1 - l{r]m-ln Qmn—l) 

T Vm—lnQ mn—1 t C m _x n C mn _x . (46) 

It is easy to check that at g = the boson variables can be integrated out in the action 
( )43|) . This may be less simple task for finite values of coupling constant go ^ 0. In the 
next section, we will apply approximations in order to eliminate completely the auxiliary 
commuting nilpotent fields from the action, and will make use of a more symmetric form 
of the integration over the bosonic fields, first over fj mn , f mn , then over rj mn , r mn . 

We would like to end this section commenting the previous exact results. We finally 
obtained a lattice field theory with action (1131) containing the same number of "fermions" 
(c, c) and "bosons" (rj, fj). Physically, this means that it is indeed possible to describe 
the system with fermionic variables for the states S = ±1 and bosonic ones for the third 
state 5 = 0. In the limit Ao — > — oo, the system is completely described in terms of 
fermions. While with A increasing to finite values, an interaction between fermions 
and bosons is added. Beyond a value Aot, fermions form bosonic pairs: in the limit 
A — > +oo, all fermions condense into bosons, leading to a purely bosonic system. In 
this interpretation, the tricritical point may be expected to be seen as a particular point 
on the critical line where the interaction is such that an additional symmetry between 
fermions and bosons appears. This might correspond to supersymmetry appearing in 
the conformal field theory describing tricritical Ising model. To our knowledge there is 
no evidence of supersymmetry derived directly from a lattice model: the exact lattice 
action ( T43l) could be a good way to see how super-symmetry may emerge from a lattice 
model. Of course all we said so far is only speculative: we are currently studying it in 
more detail, to confirm or infirm this hypothesis. 
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4. Effective action in the continuum limit 

In the Ising model case, the fermionic action on the lattice is quadratic and the 
corresponding Grassmann integral can actually be computed exactly by transformation 
into the momentum space for fermions by means of the Fourier substitution. The 
situation for the 2D BC model is less simple, as there is a non-Gaussian interaction 
part in action (PUT) , alias (l4"51) . which contains terms of order up to 8th in fermions. The 
Grassmann integral leading to the partition function can no longer be computed directly, 
as in the Ising model case, by a simple Fourier substitution. In this sense the 2D BC 
model is not integrable. However it is still possible to extract physical information by 
taking the continuous limit of the BC lattice action like (1321) . or (j4"Tj) . and analyzing it 
using tools from quantum field theory. 

4-1- Effective 2nd order fermionic field theory 

We would like to obtain an effective purely fermionic theory for the BC model up to 
order 2 in momentum k from the previous calculations, with two variables per site, to 
analyze the critical behavior of the model. In the Ising case, the critical behavior is 
given, in the continuous limit, by a massless Majorana theory that follows from two- 
variable action. In the following, we will see how to compute the mass of the BC model 
in its effective Gaussian part. The condition of the zero effective mass will give already 
the critical line in the (T, Ao) plane for BC model. For the location of a tricritical point 
on that line one needs more complicated analysis, taking into account the stability of 
the kinetic part of the action, which is in turn affected by the presence of the interaction. 
In the infrared limit, the spectrum is given by expanding the effective action, or rather 
the correspondent partial integral in Z = Yik up to the second order in the 
momentum k. The coefficient A in front of the term Afc 2 in the basic factor Zy. of Z 
is what we call the stiffness parameter of the model. It dominates all contributions 
from the kinetic part of the action. In the Ising model case, the stiffness coefficient is 
always a strictly positive coefficient. In this case, the only singularity in the spectrum 
follows from the condition of vanishing mass, resulting the Ising critical point. Here in 
the BC model, we will show that the effective stiffness coefficient can also vanish at a 
special point at the critical line in (T, Ao) plane, rendering the spectrum unstable and 
changing the nature of the singularity. This happens for large enough go, as Ao increases. 
We intend to identify the above singular point as an evidence for the appearance of a 
tricritical point, together with a segment of the first-order transition line, at the BC 
phase diagram at sufficiently strong dilution. In order to be able to perform the QFT 
analysis of the above kind, we ought to eliminate the bosonic nilpotent fields from the 
action, being interested merely in the low-momentum (small k) sector of the theory, 
and making reasonable approximations whenever necessary. 

This program also implies a more symmetric way of integration over the nilpotent fields. 
Instead of integrating over the variables r mn and f mn as in Eq. P3|) . we now proceed 
by integrating first over fj mn and f mn in Eq. (HTl) . making use of the definition of the 
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integral. This results the reduced integral with a new action: 

^ (2c COsh K^j I J"J (ic mn (lc mn (iT] mn (iT mn \c mn C mn -\- f] mn q mn TmnQmn ^]mn^~mn\ 

m,n 

x exp(5i sing + S- mt ) , (47) 
where <S>i sing is given in (|45|) . while 

«5int = gO^VmnTmn [(1 ~ 7 7 7m-ln)(l ~ 7Wl) + £ 2 C m _i n C mn _i] , (48) 

771,71 

and 

(7mn C mn C mn -\- tC mn {C m —\ n C mn — i) C mn [c mn -|- t{c m —\ n C mn — l)j , 

Omra C mn C mn -\- t('rnn(.('-rn-- in Cmn-l) [Cmn t{c m —\ n C mn _i)]c mn . (^9) 

It is also useful to note that = = 0, and q mn q mn = 0. The free-fermion Ising part 
of the action iS>i S i ng in ( 14T1) at this stage remains unchanged and is given by the standard 
expression f|45l) . The above integral (j4"Tj) includes as well the product of quadratic 
polynomial terms like c mn c mn + rj mn q mn + T mn q mn + r] mn T mn , which can not be written 
as a single exponential. However, when integrating over the remaining variables r) mn 
and T mn , it is easy to realize that these polynomial terms roughly impose the following 
substitution rules in the action iS int : 

Vnin7~mn * ^mn^mn i Vmn * QVrm ; ^rnn ^ 5mn • (^0) 

In a sense, the above rules can be considered as an operation of approximate Dirac 
delta functions on the variables 7] mn and r mn , replacing them by fermions. These rules 
of correspondence though are not unreservedly exact: when expanding the exponential 
of <Si nt into a series, the terms will appear that couple to each other to give c mn c mn but 
not q mn qmn = as is given by the above substitution rules. For example, terms such as 

(90 r ]m+ln T m+ln r y r ]mn) x (90 r ]mn+l T mn+l1 T mn ), (51) 

instead of vanishing, lead to a contribution in the effective action equal to 

do'f CmnCmnCm+lnCm+lnCmn+lCmn+l' (52) 

Therefore there are more terms in the final effective action S e s(c, c) than in the 
one resulting from the above substitution rules. However, we would like to apply 
approximations to the term with interaction, and the higher order corrections of the 
above kind can be neglected within this scheme anyhow. From an effective action that 
follows from ( 1471) . we intend to obtain the basic momentum-space factor of Z up to 
order 2 in momentum k, in order to study the stability of the free fermion spectrum. In 
the pure Ising model at criticality, the factor Z^ gives basically a (m 2 + A; 2 ) contribution 
to the partition function and free energy at small momenta, with mass m 2 = at the 
critical point [301 EH]- I n fact, there is also the stiffness coefficient A in front of k 2 in 
this term, k 2 — > Xk 2 . In the Ising case, this stiffness coefficient is non-singular at the 
critical point and can be fixed simply by its finite value at the critical temperature. 
In the BC case, however, we have a line of critical points as A varies from negative 
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to positive values. Respectively, the stiffness coefficient A = A(A ) also varies with a 
variation of the chemical potential A along the critical line. The point is that in the 
BC case the effective stiffness coefficient vanishes at some position at the critical line, 
for a sufficiently strong dilution, which may eventually be identified as the tricritical 
point of the BC model. In what follows, we apply the Hartree-Fock-Bogoliubov (HFB) 
approximating scheme [HH H?l EE] in the momentum space in order to gain a modification 
of the above Ising like behaviour provided by the presence of the interaction in the BC 
case. In essence, the HFB decouples the four-fermion interaction into few Gaussian 
terms added to the basic action. [jj] This also assumes a self-consistent calculation of 
the corrections which modify the parameters in the mass term and the kinetic part of 
the action, and eventually modify the stiffness coefficient, due to the HFB decoupling 
of the interaction. 

Among the terms that contribute in to the second order in momentum are in any case 
those coming from the kinetic part of the free-fermion quadratic piece of the action, cf. 
Eq. ()45l) . In the continuous limit, with c m _i n — ► c — d x c, c mn -\ — ► c — d y c, these terms 
are combinations of cd x c or cd x c, cd y c or cd y c. From the above rules (j5D]l . we expect 
that the effective action will contain as well quartic contributions such as ccdicdjC, 
with i,j = x,y, at the lowest order. This term is degree 4 in Grassmann variables and 
2 in derivatives. The expansion of the exponential of such terms will give corrective 
coefficients to the k 2 behaviour, and may thus change the order of the transition if the 
renormalized stiffness vanishes. We also have to consider not only the direct substitution 
of the variables with the rules given above, but also the possible correction terms like ( |52l ) 
that may contribute to the stiffness. We should also drop terms which contain a ratio of 
number of derivatives to the number of Grassmann variables higher strictly than 1/2 as 
their effect is expected to provide next-order corrections within the basic approximation 
scheme outlined above. After some algebra, these following terms contribute to the 
effective action: 

^effective Rising 90 / , CmnCmn [(1 ^Qm— In) (1 T(Zmn— l) t C m _i n C mn _iJ 



The above effective action define some lattice fermionic theory with interaction. We 
keep in mind to analyze it further on in the momentum space at low momenta, which 

| Let us remember that the interaction terms in BC model appear solely due to the presence of the 
dilute (vacancy) sites. Respectively, the strength of interaction (the coupling parameter go) increases 
with increasing rate of dilution, with variation of the chemical potential Aq. The corrections with 
go may thus appear in the mass term and the stiffness coefficients of the BC effective action within 
mean- held HFB analysis. In fact, as we shall see below, the relevant go correction to the mass at the 
Gaussian (free-fermion) level already follows when we extract the effective action from ([47)) and (150]) . 
see (|53|) . while the kinetic corrections, that at lattice level may be attributed to the correlations of the 
Ising degrees and vacancies at the same and neighbouring sites, are to be extracted self-consistently 
within HFB scheme from the residual interaction in the effective action. 



m,n 




(53) 
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corresponds to the continuum-limit interpretation of the model. 



4-2. Continuum limit 

In the continuous-limit interpretation of the above action, we replace c mn by c = c(x, y) 
and c mn by c = c(x, y), assuming as well the substitution rules like c m _i„ = c — d x c and 
Cmn-i = c — d y c. After a Fourier transformation of the fields, this corresponds to the 
low-momenta sector of the exact lattice theory around the origin k = 0. In particular, 
we put q mn — > q = cc (1 — t) + t(d x c — d y c) c and q mn — > q = c c(l — t) — tc (d x c — d y c). 
The free Ising part Rising from (l53i) gives simply 

Rising = J dxdy [(1 — 2t — t 2 ) cc — t{t + 1) cd x c + t{t + l)cd y c — tcd x c + tcd y c] . (54) 

In the above action, one can readily distinguish the mass term and the kinetic part, 
provided one assumes the QFT interpretation of the as 8 oc,ated integral^ Notice that 
the next order momentum term with product dxd y is neglected in the above action ( 1541) 
at the backbone of the first order d x , d y terms □. In the continuous limit, the last term 
of the BC action (1531) gives ccd x cd x cd y cd y c, which is 4 th order in derivatives and 6 th 
degree in Grassmann variables. The ratio of these numbers is 2/3 which is higher than 
1/2, and therefore this term can be discarded, as explained above. The term in factor 
of g m (|53l contains q m -\n and q mn -\ which need to be expanded up to the order 2 in 
derivatives, with d xx q = 2(1 — t)d x cd x c, and d yy q = 2(1 — t)d y cd y c. The effective action 
finally can be written in the continuous limit as 

S c fT = Rising + j dxdy^g cc + g cc[t(t + 2j)d y cd x c - 7(1 - t)(d x cd x c + d y cd y c)] j . (55) 

In the following, we shall use this effective action to obtain information on the phase 
diagram of the BC model. 



5. Spectrum analysis and phase diagram 

In this section we analyze the critical properties of the effective action ( l55i) and the 
low energy spectrum Zy. of Z in the momentum-space representation. In particular we 
develop a physical argument for the existence of a tricritical point on the phase diagram 
from the above fermionic action. The critical line follows already from the condition of 

% The Ising mass is easily seen from ([54| to be m lsin<r = 1 — 2t — t 2 , which must vanish at the 
critical point. Indeed, the condition of vanishing mass 1 — 2t — t 2 = gives t c = \[2 — 1, alias 
K c = J/T c = i ln(l + V2), in agreement with the exact solution of this model on a lattice. The ordered 
phase corresponds to negative mass, with t — > 1 as T — > 0. The structure of the action (|54|) rather 
implies the interpretation of the pure 2D Ising model in terms of the Majorana fermions |30, 38] [39] . 
Respectively, one may pass to the Dirac interpretation by doubling the number of fermions in the 
action. Notice also that cdiC = cdiC under the integral (|54[) since di is a skew-symmetric operator. 
+ Despite these terms with d x and d y are linear in momentum in the action, they contribute as k 2 
into the spectrum factor of Z ', while their product may only contribute at the level of next order 
corrections to Zk, as m — > 0. 
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the zero mass. At the tricritical point, we assume that the effective stiffness coefficient 
in factor Zy. also vanishes. The Hartree-Fock-Bogoliubov (HFB) approximation scheme 
will be used to count properly the effects of the interaction [HI H3, HE]. 

5.1. Phase diagram 

The BC model effective action of ( J55l) includes the free-fermion Gaussian part and the 
quartic interaction. The quadratic (Gaussian) part of the whole action is merely formed 
from iSi S i n g, but the remaining interaction term in the effective action fl55l) also includes 
quadratic term g cc, which is to be added to the Ising part of the action and will modify 
the Ising mass. The pure Ising model action is shown in (15^1) above. In the continuum 
limit, see ( T54l) . this action includes the mass term m Ising cc, with m Ising = 1 — 2t— t 2 , and 
the kinetic part. The condition for the critical point in the pure Ising case is then given 
by Ukising = [313 EH [39]. In the BC case, the presence of the Gaussian correction g cc 
will modify the mass term in the effective BC action : rn Isin — > rn BC = 1 + g — 2t — t 2 , 



which we assume to be vanishing at the critical lineJJ 
The approximations we intend to apply to tackle the remaining quartic part of the BC 
action fl55|) are of that kind that we replace, in different possible ways, the two of four 
fermions by variational parameters, or the effective binary averages, which are then 
specified self-consistently from the resulting Gaussian action. This may be viewed as a 
kind of the HFB like approximation method, which proved to be effective in systems of 
quantum interacting fermions, like BCS theory of ordinary superconductivity. This also 
implies that calculations are to be performed rather in the momentum space, but not on 
the real lattice, or its continuum real space version, and the correspondent symmetries 
are to be taken into account properly. The application of the HFB scheme also implies 
that the interaction may be not necessary weak. 

From the explicit form of the quartic part of the interaction in (|55l) . it can be seen 
that the decoupling of the quartic part of S int produces terms which only modify the 
kinetic terms in the effective action, at least, in first approximation, with calculations 
being up to order k 2 in the Zy. factor. This modification might be significant at strong 
dilution, rendering the appearance of the tricritical point and changing the nature of 
the phase transition from second to the first kind. These effects are to be discussed 
in the second part of this section. In next subsection, we consider in more detail the 
BC critical line in the (T, Ao) plane, with dimensionless temperature T and chemical 
potential A normalized by the exchange energy J. 

* The additive corrections that may contribute to the mass term from the non-Gaussian part of the 
action ([55]) are k 2 dependent and vanish as fc 2 — * 0. They may be neglected. In the effective action 
()55l) . the principal modification of the mass term due to vacancies is realized already at the Gaussian 
level by the go cc term, as it is commented above. The effect of the non-Gaussian part in (|55[) is merely 
that it produces corrections to the kinetic terms, after the HFB decoupling of the interaction. 
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5.2. Critical line 

The equation for the BC critical line we consider in this section is the one that follows 
from the condition of vanishing mass, m BC = 1 + go — 2t — t 2 = 0. In a detailed form, 
this equation reads: 



which in turn admits the explicit solution for Aq as function of T in the form: 



The inverse dependence for T as function of A can be evaluated numerically by solving 
any of the above equations, which are all equivalent to the condition of the zero mass 
in the theory with action fl55l) . This results the critical line for the BC model shown in 
Fig. 1. In the limit A — > — oo, from either of the equations (151)1) and (15"T|) . we recover the 
Ising case, with T c = 2. 269185. For finite A , as vacancies are added, we obtain a slowly 
decreasing (for moderated values of A ) function for T c = T C (A ), which terminates at 
the end-point (T c = 0, A = 2) at zero temperature, as it can be deduced from (1561) . By 
following the critical line from left to the right, at first stages, for weak dilution, from 
physical considerations and the universality argument, we expect the transition to be of 
the second kind, as in the pure 2D Ising model, while this behavior may be destructed 
for sufficiently strong dilution, as Ao increases and transfer to the positive values, where 
the correspondent term in the BC Hamiltonian already suppress the Ising states and is 
favoring vacancy states. This happens at a singular point, which we are going to identify 
from the condition of stability of the kinetic coefficient in the fermionic spectrum of the 
action (155|) . this is to be discussed in the next subsection. At the critical line, with zero 
mass, only derivative contributions remain in the action Eq. ( 1551) . These include the 
free fermion kinetic terms, and the ones presenting the residual interaction between the 
singlet level and the Ising doublet at the quartic level in fermions. The HFB decoupling 
of the interaction will modify the kinetic part of the action. 

The critical line that follows from the condition of zero mass as given by Eq. (|56l) is 
plotted in Fig. CD and compared with recent Monte Carlo simulations by da Silva et al. 
|H1 [H [9J [10] . The agreement between numerical simulations and our results is very 
good, the mass of the system (1561) being exact in that sense, at least in the transition 
region. The agreement is within 1% over the whole range of variation of A at the 
critical line T c = T C (A ), provided we use the Monte-Carlo data for T c as input and 
evaluate theoretically A from (1581) for comparison. The numerical data in the inverse 
interpretation, for T c = T C (A ) as a function of A , are given in Table 1. Note that our 
results are also compatible with exact upper bound for T C (A ) obtained by Braga et al. 
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This equation may be written as well in the form: 
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Figure 1. (color online) Comparison between critical line Eq. (|5rj)l (plain red line) 
and numerical results from Monte Carlo simulations. The black filled dots are from 
Fig. 1, da Silva et al. [9] (Wang-Landau method). The cross symbol indicates the 
tricritical point identified by the same authors. The blue diamond symbols are from 
Ref. [TO] , the magenta triangles from Ref . [7] , and the green squares from Ref. [6] (see 
also Table 1 for explicit numerical values). 



[Hj. Also notice that the value of T c (Ao) can be easily evaluated analytically at the 
point A = 0, where sinh(2/T c ) = 3/2, with the solution T c (0) = 1. 673 971 856. This is 
to be compared with the Monte-Carlo results T c = 1. 6955 ± 0.0010, T c = 1. 681(5) and 
T c = 1. 714(2) [6j [7J [JO], the agreement is again good. 

5.3. Tricritical point: Hartree-Fock-Bogoliubov analysis 

The main physical feature of the 2D BC model is the existence of a tricritical point at 
the critical line. Below this point, the phase transition goes from second order to first 
order: the tricritical point is characterized by a change in the nature of the singularity. 
This change should be seen in the BC spectrum from ( |55j) . In this section, we analyze 
the effect of the quartic terms in the action on the stability of the free fermion spectrum 
at zero mass, along the critical line go = t 2 + 2t — 1, by considering the effect of the 
interaction part of the action onto the kinetic part within the HFB like approximating 
scheme |4"B"t H71 HE]. The Ising part can be easily written in the momentum space 
representation, which we will also refer to as Fourier space, after having defined the 
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Table 1. Numerical values of the critical points (T c (Ao),Ao) in the BC model: 
comparison of different numerical simulations and equation (|56j) . Note that small 
variation of Ao causes more significant changes in T c (Aq) in the region near Ao = 2, 
as it is to be expected from (f56j) . 



following transformations: 

1 x 1 . 

c(r) = — Cfcexp(ifc.r) , c(r) = — c fc exp(-zfc.r) . (59) 

fc fc 
Using these transformations, the Ising part of the action gains block-diagonal form, 

Rising = }jt{t + l){k x - k y )(c k c h - c_ fc c_ fc ) + 2itk x c k c_ k + 2itk y c k c^ k , (60) 
fees 

where S is the set of Fourier modes that correspond to half of the Brillouin zone: if k 
is already included in S then — k is not to be included in S and vice versa (to avoid 
repetition of modes in the different sums above), so that couples of modes (fc, — fc) fill 
up the Brillouin zone exactly once. In fact, terms with fc and — fc are already combined 
together in (1601) . The mass term is dropped in flBDl) since we are on the critical line. The 
quartic term can be written in the Fourier space as 

<S in t = jp, ^2 v (k 2 ,k 4 )c kl c k2 c k3 c k4 , (61) 

ki+k'2=k:j+k4 

with the potential 

V(k 2 , fc 4 ) = -ak%k\ + a\k x 2 k% + k%k%), 

a = g t(t + 2 1 ), c/ = (7o7(l-0- (62) 

Up to now we only expressed the action in the Fourier space, or in the momentum-space 
representation, without further approximations. In order to see if the second order line 
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is stable, we use a mean-field like approximation in momentum space, similar to the 
quantum HFB method. To do so, we decompose the fourth order interacting terms 
into sums of quadratic terms with coefficients to be determined self-consistently. These 
coefficients are actually two-point correlation functions for fermions in the momentum 
space. The interaction can be decoupled in different ways. For example, considering 
the terms contributing to the Ising action, we may take account of the averages (ckCk) , 
(c_fcC_fc), (c fc c_fc) and (cfcC_fc). There are also three different ways to decouple the 
interacting term, since can be paired with either of Ck 2 , Ck 3 , or Ck 4 - For example, 

CfcjC^ = (c fcl Cfc 2 ) + (c fcl Cfc 2 - (c fcl Cfc 2 )) = (c fcl c fc2 ) + 5 ClC2J (63) 

where 5 ClC2 is assumed to be a small fluctuation. In this case, from Eq. (160]) . the average 
is non zero only for k\ = — &2 = k or — k, with k e S. We can pair the other terms 
by writing the action in the g$ = 3 different possible ways that are compatible with the 
symmetries of Eq. (1601 . and by using the fermionic rules, we write: 

Sint = y^— y( k 2,k A ) ((c kl c k2 ) + 5 ClC2 ){(c k:j c ki ) + 5 5:i c 4 ) (64) 

ys fc 1 +fc 2 =fc 3 +fc 4 

-((CfciCfc 3 ) + <5 Cl c 3 )((Cfc 2 C fc4 ) + 5 C2 c 4 ) + ((c fcl C fc4 ) + 5 Cl c 4 )((Cfc 2 C fc3 ) + 5 C253 ) . 

The next step is to discard terms that are proportional to the squares of fluctuations 5 2 , 
and keep the others. After some algebra, we obtain the mean-field quadratic operator 
for the interaction term as follows: 

«£nt = T2 — ^ c kC-k(c k 'C^ k ')V(k, k') +4c fc c_ fc (c fc /c_fc/)y(fc / , k) 
@ s fc.fc'es 

+c k c k {{ck'Ck')v(k,k') + (c_ fc /c_ fe /)?;(fc, -fc')) 

+c^ k C-k{{ck>Ck')v{-k, k') + (c_fc/c_ fc /)i;(-fc, —k')\ (65) 
where we have defined the potential 

v{k, k') = -V(k, k) - V(k', k') + V{k, k') + V(k', k). (66) 
In the above expressions, there are three different kinds of quantities, that contribute to 
the action, associated with the sums like 2^fc c fcCfc, J2 k Ck ^ k ^^ or J2k Ck ^ k ^j^ with 
i, j = x,y. The first term gives a contribution to the total mass, the second one 
corresponds to current operators, and the third one can be thought as a dispersion 
energy tensor. Considering the symmetries of the Ising part, and the fact that the 
action must be invariant by a dilation factor at criticality, we may only take into account 
the current operators. Respectively, we can drop the first two terms in the potential 
v(k, k') defined in Eq. ( |66i) . We define therefore the following unknown parameters, for 
the diagonal and nondiagonal couplings of fermions (i = x, y) : 

k = 777^y~l(( c fcCfc) - (c-kC-k))h, 
kes 

Ui = Jjl ^2( c kC-k)ki, Ui = ^(c fc c_fc)/cj. (67) 
fees kes 
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From the previous discussion, we can drop the first two terms in the potential v(k, k') 
defined in Eq. fl66|) . since we already assume that only currents are kept as parameters 
along the critical line. In this case, it is easy to rewrite, from the property v(k, k') = 
—v(—k, k') = —v(k, —k'), the effective mean field action of (1631) as: 



S int = — ^ Aic k c- k [(au y - a'u x )k x - a'u y k y ] + Aic k C- k [-a'u x k x + (au x - a'u y )k y ] 
gs t~s 

+2i(c k c k - c^ k c- k )[(at y - 2a't x )k x + (at x - 2a't y )k y }. (68) 

We make then further assumption that, by symmetry invariance in the momentum 
space, there exists a solution satisfying u y = u x , u x = u y and t x = —t y , so that: 

•Sint = — y^4zc fc c_fc[(att a; - a'u y )k x - a'u x k y ] + Aic k c^ k [-a'u x k x + (au x - a'u y )k y ] 



9s 



kes 



-2i(c fc c fc - c_ k c_ k )(k x - k y )(a + 2a')t x . 

The total effective action (with zero mass) can finally be written as 

2 



(69) 



E 



kes 
9s 



t(t + l) 



9s 



{a + 2a')t x 



(k x ky)(c k c k c^ k C- 



+i ~ I f^T* + ( au x ~ a ' u y)) k x - a'u x k y c k .c- k 
9s 

+i— -a'u x k x + (^-t + (au x - a'uyfj k y c k C- k , 
or in a more compact form as 

S cS = 22 ic (kx ~ K)( c kCk - C-fcC_fc) + 2i(ak x - bk y )c k C- k 
kes 

+2i(-bk x + ak y )c k c_ k , 
with the following coefficients 



(70) 



a = t + 2- 



au T — olu 



2a 



c = t(t + 1) - 2t 



a + 2a' 



(71) 



(72) 



9s 9s 9s 

The partition function can then be written as a product over the Fourier modes 

z = U k es Z ^ with 

Z k = k 2 [A + Bsm29 k ], (73) 
9k being the angle of the vector k, and 

A = c 2 - Aab, B = -c 2 + 2{a 2 + b 2 ). (74) 

We assume that \A\ is larger than \B\ on the second order critical line, until a singular 
point is reached, where eventually A 2 = B 2 . Indeed, the expression f !73p is valid only if 
the elements A + Bsin26k are all strictly positive, which is the case only if A 2 > B 2 . 
This will be checked using numerical analysis. Beyond this point, the effective action 
is unstable and has to be modified to incorporate further corrections. In a bosonic 
$ 6 Ginzburg-Landau theory describing a first order transition, the tricritical point is 
usually defined as the point where both coefficients of $ 2 and $ 4 terms vanish 
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By analogy, in the present fermionic theory, it is tempting to associate the above singular 
point with the effective tri critical point. 

The parameters t x , u x and u y are to be determined self-consistently from the definitions 
Eqs. ( 1671) . In the continuous limit, these reduce to 



Numerically, we proceed the following way. Starting from T slightly below T c (— oo), 
we solve the consistency equations for t x , u x and u y , with the value of A given by the 
critical line ( l56l) at a given temperature. The solutions are then plug into the coefficients 
A(T) and B(T), and we plot A(T) 2 — B(T) 2 as a function of T, as is shown in figure [2] 
We repeat the process by decreasing the temperature until we reach the point where 
this quantity vanishes. 

By doing so we find a singular point approximately located at (T t *,AQ t ) ~ 
(0.42158,1.9926). This is close to the tricritical point T t given by Monte Carlo 
simulations: (T t ,A , t ) ^ (0.610,1.9655) [9J, and (T t , A , t ) ^ (0.609(3), 1.966(2)) [TO]. 
If we assume that T* represents the tricritical point, the mean- field like treatment of 
the underlying field theory underestimates the fluctuations, rendering the second order 
critical line more stable at lower temperatures, as compared to Monte-Carlo results, 
as we approach (T c = 0, A = 2) along the critical line. Stronger fluctuations can be 
simulated by lowering the value of gs, which increases (lowers) the value of T* (AqJ, 
respectively. Instead of gs = 3, taking gs = 2.5, for example, leads to a T* ~ 0.48, 
closer to the Monte Carlo results. This can be achieved precisely by incorporating more 
diagrams in the computation of the effective free energy [47]. Also, due to the fact that 
we are in a region near (T c = 0, Ao = 2), where the change in temperature is large 
compared to the change of A (the slope is vertical at this point as is seen in figured]), 
it is more difficult to obtain a precise value of T t * within a mean- field treatment. 
It is important that the BC fermionic action (!55j) finally predicts the existence of a special 
(tricritical) point at the critical line somewhere close (in Ao) to termination point of 
that line at (T c = 0, Ao = 2). The tricritical point is defined, within this interpretation, 
as the point of the destruction, or loss of stability, in the effective fermionic spectrum 
of the action due to the modifications introduced into the kinetic part by a sufficiently 
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After computing the trigonometric integrals, we obtain the relations 
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Figure 2. (color online) Stiffness of the spectrum: solution of HFB self-consistent 
Eqs. (T?6"]) for the coefficient A(T) 2 — B(T) 2 as function of T. The temperature where 
A(T) 2 - B(T) 2 = gives the location of the singular point T*. 



strong dilution of a system by the vacancy states, which corresponds to large enough 
coupling constant go, as h was commented above. |]J. 

6. Conclusions 

In this paper, we have considered the physics of the BC model as a fermionic field 
theory. Using Grassmann algebra, we have shown that the model can be transformed 
into quantum field theoretical language in terms of fermions alias Grassmann variables. 
This fermionic theory for BC model is described by an exact fermionic action with an 
interaction on a discrete lattice. This action can be reduced, after some transformations, 
in the continuum limit and low energy sector, to an effective continuum field theory 
which includes a modified Ising action, which is quadratic in fermions, and a quartic 
interaction. From there we have extracted the exact mass of the model and analyzed the 
effect of the quartic term on the stability of the free fermion spectrum in the kinetic part. 
The condition of the zero BC mass gives the critical line of phase transition points in the 
(T, A ) plane, which is found to be in a very good agreement with the results of Monte- 
Carlo simulations over the whole range of variation of concentration of the non-magnetic 

j It may be also noted that the Monte-Carlo values for (Tt,Ao,t) seemingly lie practically on the 
theoretical curve for the critical line ([56]) - ([58]) . For instance, taking as input value T t ~ 0.609(3) [ID] , 
from ([58)) we find Ao.t — 1.952, which is sufficiently close to the M-C value Ao.t — 1.966(2) from this 
set [10], the deviation being probably less than 1%. 
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sites governed by A . The location of the tricritical point needs additional analysis of 
the excitation spectrum of integral factors of Z around the origin in the momentum 
space. In particular, the stiffness of the excitation spectrum (the coefficient in front of 
k 2 term in factors Zy. as we expand the dispersion relation for Z in momentum variables) 
vanishes at a singular point T*, which we assume to be identified as the tricritical point 
T t . A Hartree-Fock-Bogoliubov analysis gives an approximate location for this point 
on the phase diagram (critical line) which can be compared to the numerical results 
of Monte Carlo simulations. The more precise location of the instability point could 
be achieved by taking into account more diagrams contributing to the effective free 
energy. In any case, we have shown the existence of a singular point at the critical line 
by studying the stability of the kinetic spectrum of the action at this line, where the 
nature of the transition is to be changed due to strong dilution. The main result of this 
paper is the possibility to study precisely first-order transition driven systems from a 
fermionic point of view using Grassmann algebra. The method we have applied may 
be useful as well for other systems where effective field theory is presented by an action 
similar to that of Eq. fl55|) . In essence, this is a one of the simplest form of an action 
with 4-fermion interaction that can be written out from a unique pair of Grassmann 
variables at each point of the real space in two dimensions. Application of the same 
method to other extensions of the BC Hamiltonian, such as the Blume-Emery-Griffiths 
model [3] , is also possible. Finally, at intermediate stages, a partial bosonization of the 
system leads to a mixed representation of the model not only in term of fermions but 
also in term of hard core bosons, as written explicitly in the lattice action of Eq. ( 1431) . 
The representations of this kind could be useful also to look for a possible interpretation 
of the tricritical point in the BC model as a special point in the phase diagram where 
an additional hidden symmetry between fermions and bosons may appear. 
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